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We comment on the relation of our previous work on the classical gluodynamics of high energy nu- 
clear collisions to recent work by Lappi (hep-ph/0303076 ). While our results for the non-perturbative 
number liberation coefficient agree, those for the energy disagree by a factor of 2. This discrepancy 
can be traced to an overall normalization error in our non-perturbative formula for the energy. When 
corrected for, all previous results are in excellent agreement with those of Lappi. The implications 
of the results of these two independent computations for RHIC phenomenology are noted. 

PACS numbers: 24.85.+p,25.75.-q,12.38.Mh 



I. INTRODUCTION 

In a series of recent papers, we computed numeri- 
cally the classical gluodynamics of the early stages of 
very high energy nuclear collisions. These computations 
were first performed for an SU(2) Yang-Mills gauge the- 
ory PL Gill an d later for the SU(3) Yang-Mills gauge 
theory 0,1a, @- In a very recent paper, Lappi has inde- 
pendently performed the same computation for the SU(3) 
gauge theory 0. We will show that the apparent dis- 
crepancies in the two computations can be understood 
as arising from an incorrect overall normalization in the 
non-perturbative formula for the energy in our previous 
work. Since the problem arises from a normalization er- 
ror, all of our previous results are still useful if interpreted 
appropriately. The purpose of this note is primarily to 
clarify the sources of discrepancy in the two computa- 
tions and to present corrected results wherever necessary. 
Now that two independent computations can be shown 
to agree, we will show that they significantly constrain 
phenomenological interpretations of the RHIC data. 

This note is organized as follows. In Section II, we 
discuss the different conventions for the small x classical 
effective Lagrangean density used in the literature and 
isolate the source of the normalization error in our pre- 
vious computations. A reader uninterested in the details 
can skip right ahead to Section III where we discuss the 
corrected results for various quantities. The implications 
of these results for RHIC phenomenology are discussed 
in Section IV. 



II. NORMALIZATION OF CLASSICAL 
EFFECTIVE LAGRANGEAN 

There are two conventions for the classical effective 
Lagrangean density 0] that are followed in the literature. 
This is occasionally confusing so we shall discuss these 



explicitly here. The first of these is 

L I [A,p) = - A F a F a -p a A a +i^ r . (1) 

For convenience, all Lorentz indices and space-time in- 
tegrals have been suppressed here. For example, F a de- 
notes the field strength tensor. In this form of the ef- 
fective Lagrangean density, the coupling constant g is 
absorbed in the definition of the charge density p a . The 
color charge squared per unit area p\ contains a factor 
g 2 . The correlator of color charges is then 

(p a (x ± )p a (y ± ))=p 2 A S ab 5^(x ± ~y ± ). (2) 

This form of the effective Lagrangean density is used in 
several works - for instance, Refs. [flLUcj. 

The other form of the effective Lagrangean density 
found in the literature, for instance in Refs. |8l Ull Il2] ] is 

L II [A,p\ = -F a F a -gp a A a + i^-. (3) 

Here g is not absorbed in the definition of the charge 
density p a but appears explicitly in the second term of 
the classical effective Lagrangean density. The relation 
between the color charge squared values in the two La- 
grangean densities is therefore 

2 22 
Ma = 9 • 

The form of the effective Lagrangean density in Eq. © 
is the form used in our papers [E 0. S, S B Q j and the 
scale we use, A 2 , is defined to be A 2 . = ,g 4 /i 2 v . Eq. (1) of 
Ref. has the form 

(p a (x x )p a (y ± )) =g 2 p, 2 N 6 ab SW(x ± -y ± ). (4) 

Thus from Eqs. J3J and 10} , it is clear that the form of 
the effective Lagrangean density used by Lappi is Eq. JQ), 
albeit with the color charge squared defined to be iden- 
tical to ours. 
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In Eq. and @, the field strength tensor has the 
canonical form 

F; v = d^Al - d v A% + gf ahc A\Al . (5) 

In our numerical work - see for instance Ref. Q - the 
field strength tensor is defined as 

Fl» = d X ~ d v Al + 2f abc A b X ■ (6) 
Performing the field re-definition A° = ^A®, we obtain 



where f E was determined by an extrapolation of the lat- 
tice results to the continuum limit. Then from Eq. (|lf)|> 
and JEJl, we obtain 



4 -f t 2— \3 



(14) 



However, ^ w = /ijv/2. Performing this substitution, we 
obtain 



1 * A3 
= -r f E A* , 

g 2 



(15) 



(7) 



The "canonical" effective Lagrangean density in Eq. 
can then be re-written as 



^ii 



\4 2 

t - „ . 
2/4 



(8) 



With a further re-definition p a = g p a /2, we obtain 



4 / 1 -— a— — u ^— r-a 



(9) 



where 7% = MJv/2. With the field re-definition, A a — > 
A a /g, both Lj and Lu, as defined in Eqs. and @ 
respectively, can be expressed in the same form as the 
right hand side of Eq. J^J , albeit without the factor of 4 
in front of the first two terms, and with /ijv instead of 
~p N . Since the continuum limit of the lattice Hamiltonian 
in Ref. corresponds to Lj with A a — > ^4 a /<7 and p a = 
/3 a /g, there is this overall difference of a factor of 4 in the 
normalization between our lattice Hamiltonian and that 
of Lappi plus the requirement that the charge densities 
are related by ~p N = /ijv/2. 

Now consider the expression for the energy per unit 
rapidity per unit area, 



et = i (E a E a + B a B a ) , 
again suppressing Lorentz indices. From Eq. 0, 



which gives 



£T — — £T . 

g 2 



In Refs. and [H, we computed 



(10) 



(11) 



(12) 



(13) 



where /s = /e/2. Thus the physical value of the energy 
liberation coefficient is one half that computed previ- 
ously. It accounts precisely for the result found by Lappi. 

The number liberation coefficient remains the same af- 
ter the field re-definitions. This is because the number 
per unit rapidity per unit area is defined to be 



(16) 



Since ~p N = /itjv/2, this replacement here exactly can- 
cels the extra factor of 4 appearing in the normalization. 
Hence f N = /at - our previous result is in agreement 
with the result of Lappi. 

The gluon number distributions in the two papers coin- 
cide if we replace A s — > A s /2. The resulting comparison 
of the result of Lappi to our prior result in Ref. |4| is 
excellent and is shown in Fig. ^ 



III. CORRECTED RESULTS 

Now that we understand the origin of the apparent 
discrepancies in the numerical results of Ref. Q and our 
prior results, we will discuss the corrected results here. 
As discussed previously 0, all dimensional quantities in 
the classical effective theory can be expressed in terms 
of the appropriate power of A s times a non-perturbative 
function of A S R. This is because A s and R are the only 
dimensional scales in the problem. A s is, in principle, 
a function of the energy, the centrality and the atomic 
number. For RHIC energies, one can broadly estimate 
that A s ~ 1-2 GeV (more on this in the next section). 
The radius of a gold nucleus is approximately 6.5 fm, so 
A S R lies in the range 32.5 - 65 for Au-Au collisions at 
RHIC. We obtain, for the transverse energy, 



1 dE T 



irR 2 di] 



f E (A s R) 3 

g 2 



(17) 



Previously, in Ref. Q, we had / B (25) = 0.537 and 
/ B (83.7) = 0.497. Following Eq. (JT5JI, this should be 
corrected to read /is (50) = 0.27 and /b (167.4) = 0.25. 
The corresponding expression for the gluon number is 



1 dN 
^Ri ~d~r) 



f N (A s R) A2 



(18) 
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with /jv(50) = 0.315 and /W(167.4) = 0.3. The variation 
in the number with A S R is very small, on the order of 
10%. These results are in good agreement with those of 
Lappi. 

As we have explained in detail in our earlier papers, 
there are a variety of ways of defining the particle number 
for a generic field <f> and its conjugate momentum 7r. Let 
us compare the results for the following three definitions 
(here and in the following kx is the transverse momentum 
of the field): 



n(k T ) 



^(MMMI 2 } 

(|7T(fc T )| 2 )Mfc T ). 



(19) 
(20) 
(21) 



Here co(kT) is the eigenfrequency for a free massless plane 
wave of the wave number fey on a square lattice. In our 
work, as in Lappi's, these prescriptions are applied to 
field configurations in the Coulomb gauge. In Fig.n° ur 
transverse momentum distributions for gluons are com- 
pared with Lappi's result. KNV I (circles) in Fig.^is the 
same distribution as in our previous papers, but with the 
correct value of A s used to set the horizontal scale. KNV 
I corresponds to the gluon number distribution from the 
definition Eq. I|19|l .and kr dependence is extracted only 
for the transverse momenta along the principal lattice 
directions. The solid line is from Lappi's numerical re- 
sult, which is computed according to Eq. (|21JI using the 
entire Brillouin zone. KNV II is obtained by the same 
definition as Lappi's and it is consistent. The deviation 
at large kx is considerd to be a consequence of how we 
define kr on the lattice. It should be noted that, with the 
lattice spacings used, the ultraviolet portion of the spec- 
trum, kx/A s > 3, cannot reliably reproduce continuum 
physics. 

The initial transverse energy per particle changes (be- 
cause Je does) and is now 



~N 



0.9A S 



(22) 



This value is nearly a factor of 2 lower than previ- 
ously |2g. As we will discuss in the following section, this 
lower value of the transverse energy per particle makes 
the connection to RHIC phenomenology simpler to inter- 
pret. 

The "formation time" Tj (defined as = 
1/^(A S R)/A S ) of gluons extracted in Ref. 0] (it is 
the same for SU(2) and SU(3) as confirmed in Ref. 0) 
is unaffected, albeit the A S R at which it is evaluated 
should be a factor of 2 larger than stated. One has 
approximately 7 = 0.3 in the range of interest. The 
initial energy density is then 



1 dE T 



7rR 2 dr) 



0.08 



A 



(23) 



The number distributions, with coefficients appropri- 
ately scaled with A s , have the same form as previously: 
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FIG. 1: Comparison of gluon transverse momentum distribu- 
tions per unit area as a function of kr/As. KNV I (circles): 
the number defined as in Eq. 119L with kr taken to mean 
the lattice wave number along one of the principal directions. 
KNV II (squares) and Lappi (solid line): the number defined 
as in Eq. 12 U . with an average over the entire Brillouin zone 
and with kr taken to mean the frequency w(fcr). 



1 



dN 



irR 2 drjd?kT 



r/n(fcr/A s ), 



(24) 



where f n (k T /A s ) is 



01 



In 



exp 



( y/k» + m 2 /T eff ) - l] 1 (k T /A s < 1.5) 



a 2 A4log(47rA;T/A s )fc T 4 



= 0.0358A. S , T off = 



{k T /A s > 1.5) 
(25) 

0.465A,, and 



with ai = 0.137, 
a 2 = 0.0087. 

With this expression for the number distributions, we 
can now compute the occupation number for gluons af- 
ter the collision. The relevant quantity for the validity of 
the three dimensional classical approximation is the three 
dimensional occupation number [27J which can be esti- 
mated from the two dimensional boost-invariant number 
distribution computed in our simulation by the following 
relation 1131 : 



glue 



(27T) 3 



1 



dN 



2 (N 2 - 1) irR 2 dr]d 2 k T 



(26) 



Substituting our result for the number distribution from 
Eq. I|19(l here, one can compute the occupation number 
of gluons at the early stages of the collision. The results 
are shown in Fig. [21 for A S R = 66. Strictly speaking, 
the classical description is valid for / > 1. It is not 
known however what the value of / is at which the clas- 
sical description breaks down completely and quantum 
corrections are dominant. 
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FIG. 2: Gluon occupation number denned in Eq . 1261 . com- 
puted using Eq. I|19^ . 



system [l^. Moreover, if the system were to approach 
thermalization due to re-scattering, this process would 
again cause the initial transverse energy to be transferred 
into longitudinal energy . Even if there were no pdV 
work done and no re-scattering, it is difficult to envi- 
sion a scenario where the final transverse energy is larger 
than the initial one. The second constraint is also very 
plausible. For a thermal system, entropy is simply pro- 
portional to the multiplicity and the final entropy must 
be equal to or greater than the initial one. Even in the 
other extreme scenario, independent fragmentation with 
no re-scattering, parton-hadron duality suggests that the 
parton multiplicity prior to hadronization is equal to the 
multiplicity of hadrons after ^3] • 

Combining the first condition in Eq. (|28J) on the trans- 
verse energy with Eq. I|17|l and the empirical result for 
the RH1C data, we find that 



A* HIC > 1.3 GeV. 



(29) 



All of our other results can be understood by appro- 
priate re-scaling. For instance, the pt distribution of the 
azimuthal asymmetry, V2(pt), plotted versus pt/A s is the 
same after replacing A s — > A s /2 on the x-axis of the 
plot p|. The peak in V2{pt) is therefore at even lower 
momenta then previously - at pt ~ A s /8. 



Here we have assumed that the average transverse area 
for the range of impact parameters corresponding to the 
most central collisions is irR 2 = 130 fm 2 and g = 2 j2^| 
Combining the second condition in Eq. (|28|l on the par- 
ticle multiplicity (at 77 = 0) with Eq. I|18|) . we find 



A* HIC <2GeV. 



(30) 



IV. NUMERICAL GLUODYNAMICS AND 
RHIC PHENOMENOLOGY 



Also, from Eq. J22J), and Eqs. we obtain the 

following constraint for the initial transverse energy per 
particle: 



Now that we have shown that the results of two inde- 
pendent computations converge, it is useful to consider 
its implications for RHIC phenomenology. The results 
obtained in this classical approach are valid only at the 
very early stages of the nuclear collision and the connec- 
tion to the final state distributions of measured hadrons 
might be considered remote. We will argue that, on the 
contrary, very general arguments can be applied to the 
RHIC data which will help constrain the properties of the 
initial state. In addition, they also significantly limit the 
freedom of final state interaction models which attempt 
to describe thermalization and the subsequent evolution 
of a Quark Gluon Plasma (QGP) in very high energy 
heavy ion collisions. 

We know that the hadronic multiplicity at rj « for 
central Au- Au collisions with center of mass energy y/s = 
130 GeV/nucleon at RHIC is ~ 1000 while the hadronic 
transverse energy is ~ 500 GeV Q. One expects, on 
general grounds, the following relations, respectively, to 
hold for the transverse energy and the multiplicity: 



Ex 



RHIC 



E 



eluc 



> E 



N glue < N 



T 

hadrons 



(27) 
(28) 



The inequality for the transverse energy follows from 
the expectation that some of the initial transverse en- 
ergy is converted into work due to the expansion of the 



1.14 GeV < — - < 1.76 GeV, 

iV initial 

and, for the initial energy density, (from Eq. (|23|l '). 



GeV DTTjn GeV 
7.1 — < S < 40 



fm 3 



fm 3 



(31) 



(32) 



Thus the experimental data, the non-perturbative for- 
mulae in Eqs. (|17|) - H18|) obtained from the classical nu- 
merical gluodunamics, and very general constraints can 
all be put together to significantly restrict the allowable 
range of A™ 10 . These in turn restrict the initial trans- 
verse energy per particle and the initial energy density. 
How these quantities evolve tell us a great deal about 
the space-time evolution of partonic matter in a heavy 
ion collision. 

In our analysis, A s is an external parameter - A 2 = 
g 4 /i 2 , where [j, 2 , as discussed previously, has the physical 
interpretation of being the average color charge squared 
per unit area per unit rapidity of color sources at higher 
values of x than those of interest. This scale is not quite 
the same as the gluon saturation scale Q 2 , which de- 
notes the scale at which the gluon distribution changes 
qualitatively due to saturation effects 9j. The relation 
between the two can be obtained by comparing the ana- 
lytical expressions for the gluon distributions in the two 
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approaches 0,0,0. One obtains 

For the values of A s in the range of interest, Q s ps A s . 
Extrapolating the Golec-Biernat-Wusthoff parametriza- 
tion to RHIC (with the appropriate A 1 / 3 dependence), 
we find A 2 sa Q 2 S = 2 GeV 2 . This value for Q s is also ob- 
tained by Kharzeev and Nardi |20j . A s thus obtained lies 
within the range required by Eqs. I|29|l and Eq. H30(l . Nev- 
ertheless, there is an ambiguity in Eq. I|33l) . For smaller- 
values of x, one expects the infrared scale in this equation 
to be of order Q s and not Aqcd [llluSf: the relation be- 
tween Q s and A s will then be modified. Secondly, even at 
the classical level, the requirement that color neutrality 
be maintained |23j will alter the relation between A s and 
Q s . Indeed, in this case, one can directly relate the two 
scales with an appropriate definition of Q s . This point 
will be discussed further at a later date. Understand- 
ing this relation is important because the constraint on 
A s discussed here will translate into a constraint on the 
gluon content in the nuclei at RHIC energies. 

The bounds in Eq. I|31|) and Eq. (|32[) are very impor- 
tant for understanding the "late" stage dynamics of high 
energy heavy ion collisions, namely, when the classical 
approach breaks down as it must. One possibility is that 
the system undergoes hydrodynamic expansion [lfij . At 
the end of the classical stage |2^, the system is com- 
pletely out of equilibrium. This is because the momen- 
tum distribution is very anisotropic: p t ~ A s and p z ~ 0. 
If the system is to thermalize, there must be particle 
production after the classical stage. Particle number con- 
serving 2 — >■ 2 interactions are too inefficient to thermal- 
ize the system within a reasonable time scale ^| |24j . If 
A s is 2 GeV, particles cannot be produced subsequent 
to the classical stage because this value saturates our 
bound for the particle number. If A s is 1.3 GeV (at 
the lower end of our allowed range), one can certainly in- 
crease the particle number by approximately a factor of 2 



after the classical stage. The transverse energy per parti- 
cle (see Eq. JSD) will then also decrease by a commensu- 
rate amount and approach the experimental value. How- 
ever, one cannot then have any hydrodynamic expansion 
(and therefore thermalization) because for A s = 1.3, the 
transverse energy was already as low as it could be: hy- 
drodynamic expansion would lower it still further. 

The problem with this "free streaming" scenario is that 
one does not have a satisfactory mechanism to explain 
the momentum dependence of the azimuthal anisotropy 
(especially at low pt) observed at RHIC. The classical nu- 
merical simulations of V2 do not generate enough V2 0- 
There is an interesting attempt to explain the RHIC data 
as resulting from "non-flow" correlations but at low p t 
(pt < Qs) its energy dependence is different from the 
trend in the measured data from RHIC at \fs = 130 
GeV/nucleon and ^/s = 200 GeV/nucleon j2£j. A con- 
sistent picture of all the global features of the RHIC data 
might still be feasible within the straitjacket provided by 
the classical simulations and the data but there is limited 
room for models of final state interactions to maneuver. 



Acknowledgments 

We would like to thank Rolf Baier, Keijo Kajantie, 
Dima Kharzeev, Tuomas Lappi, Larry McLerran and Al 
Mueller for very useful discussions and Tuomas Lappi and 
Keijo Kajantie for correspondence on their work. R. V.'s 
research was supported by DOE Contract No. DE- 
AC02-98CH10886 and the RIKEN-BNL Research Cen- 
ter. Y. N.'s research is supported by the DOE under Con- 
tracts DE-FG03-93ER40792. A. K. and R. V. acknowl- 
edge support from the Portuguese FCT, under grants 
POCTI/FNU/49529/2002 and CERN/FIS/43717/2001 
and support in part of NSF Grant No. PHY99-07949. 
A. K. and R. V. also wish to thank the INT at the Uni- 
versity of Washington for partial support during the com- 
pletion of this work. 



[1] A. Krasnitz and R. Venugopalan, hep-ph/9706329 

hep-ph/9808332 Nucl. Phys. B557 237 (1999). 
[2] A. Krasnitz and R. Venugopalan, Phys. Rev. Lett. 84 

(2000) 4309. 

[3] A. Krasnitz and R. Venugopalan, Phys. Rev. Lett. 86 

(2001) 1717. 

[4] A. Krasnitz, Y. Nara and R. Venugopalan, Phys. Rev. 

Lett. 87, 192302 (2001). 
[5] A. Krasnitz, Y. Nara and R. Venugopalan, Phys. Lett. B 

554, 21 (2003); |arXiv:hep-ph /0209341 
[6] A. Krasnitz, Y. Nara and R. Venugopalan, Nucl. Phys. 

A 717, 268 (2003). 

[7] T. Lappi, |arXiv:hep- ph/0303076"l 

[8] L. McLerran and R. Venugopalan, Phys. Rev. D49 2233 
(1994); D49 3352 (1994); D50 2225 (1994). 



[9] Y. V. Kovchegov and A. H. Mueller, Nucl. Phys. B 529, 
451 (1998). 

[10] E. Iancu, A. Leonidov and L. McLerran, 
arXiv:hep-ph/0202270 E. Iancu and R. Venugopalan, 
|arXiv:ne p-ph/0303204 

[11] J. Jalilian-Marian, A. Kovner, A. Leonidov, and H. 
Weigert, Nucl. Phys. B504 415 (1997). 

[12] M. Gyulassy and L. McLerran, Phys. Rev. C56 (1997) 
2219. 

[13] A. H. Mueller, |arXiv:hep -ph/0208278 

[14] . Adler et al. [STAR Collaboration], Phys. Rev. Lett. 87, 
112303 (2001); K. Adcox et al. [PHENIX Collaboration], 
Phys. Rev. Lett. 87, 052301 (2001); I. G. Bearden et 
al. [BRAHMS Collaborations], Phys. Lett. B 523, 227 
(2001); B. B. Back et al. [PHOBOS Collaboration], Phys. 



6 



Rev. C 65, 061901 (2002). 
[15] K. J. Eskola, P. V. Ruuskanen, S. S. Rasanen and 

K. Tuominen, Nucl. Phys. A 696, 715 (2001). 
[16] A. H. Mueller, Nucl. Phys. B572 (2000) 227; J. Bjoraker 

and R. Venugopalan, Phys. Rev. C63 024609 (2001). 
[17] Ya. I. Azimov, Yu. L. Dokshitzer, V. A. Khoze, and S. I. 

Troyan, Z. Phys. C27:65, (1985). 
[18] J. Jalilian-Marian, A. Kovner, L. McLerran and H. 

Weigert, Phys. Rev. D55 5414 (1997). 
[19] Y. V. Kovchegov, Phys. Rev. D 54, 5463 (1996). 
[20] D. Kharzeev and M. Nardi, Phys. Lett. B 507, 121 

(2001); D. Kharzeev and E. Levin, Phys. Lett. B 

523, 79 (2001); D. Kharzeev, E. Levin and M. Nardi, 

|arXiv:hep-ph/01113T5] 
[21] E. lancu, K. Itakura and L. McLerran, Nucl. Phys. A 

708, 327 (2002); 

[22] A. H. Mueller, arXiv:hep-ph/0206216 

[23] C. S. Lam and G. Mahlon, Phys. Rev. D 61, 014005 



(2000); ibid., 62, 114023 (2000). 
[24] R. Baier, A. H. Mueller, D. Schiff and D. T. Son, Phys. 

Lett. B502 51 (2001); B539, 46 (2002). 
[25] Y. V. Kovchegov and K. L. Tuchin, Nucl. Phys. A 708, 

413 (2002). 

[26] It is not exactly a factor of 2 because the value of A 3 R 
at which /e is evaluated changes as well. 

[27] We thank Al Mueller for a discussion on this point. 

[28] There is clearly an ambiguity in this bound due to our 
choice of the strong coupling constant g, which cannot be 
resolved at this order of the computation. If we assume 
that a 3 = a s (A s ), the one loop value of a s (from the 
bound we will set) is close to this estimate. 

[29] It is estimated in Ref. [24| that the average occupa- 
tion number /, for a s <C 1, is less than unity at r > 



